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VARIATIONAL CALCULATIONS FOR RESONANCE OSCILLATIONS 
OF INHOMOGENEOUS PLASMAS 

b y 

Y-K. M . Peng and F. W. Crawford 
Institute for Plasma Research 
Stanford University 
Stanford, California 94305 

ABSTRACT 

In this paper, the electrostatic resonance properties of an inhomo- 
geneous plasma column are treated by application of the Rayleigh-Ritz 
method . In contrast to Parker, Nickel, and Gould (1964), who carried 
out an exact computation, we have used a description of the rf equation 
of motion and pressure term that allows us to express the system of 
equations in Euler-Lagrange form. The Rayleigh-Ritz procedure is then 
applied to the corresponding Lagrangian to obtain approximate resonance 
frequencies and eigenfunctions . An appropriate set of trial coordinate 
functions is defined, which leads to frequency and eigenfunction esti- 
mates in excellent agreement with the work of Parker, et al. (1964). 



1 . INTRODUCTION 


This paper is concerned with the use of the Rayleigh-Ritz procedure 
to estimate the electron resonance frequencies of a warm inhomogeneous 
plasma column. This procedure has been extensively applied to single self- 
adjoint equations with great success (Mikhlin, 1964). For a system of 
equations, however, theoretical extensions have been noted only for the 
case of elliptic equations (Mikhlin, 1965), For the equations to be 
used here, which are not elliptic, it will be shown that accurate resonance 
frequencies can be predicted, provided that a certain set of coordinate 
functions is defined. 

Previous theoretical treatments of the electron resonance problem 
give predictions which agree well with experimental data. However, these 
approaches have encountered difficulties stemming from the inhomogeneous 
electron density profile. For example, Parker, Nickel, and Gould (1964) 
solved numerically an appropriate fourth-order differential equation for 
the rf potential in a cylindrical positive column. Because of the exponential 
nature of the solutions in the cutoff region, their calculations were 
limited by the condition, a 2 /\ 2 < 4500, where a is the column 
radius, (= e^rcT/n^e 2 ) is the mean-squared Debye length, and n^ 

is the mean electron density in the column. Baldwin (1969) used a 
kinetic model in the low- temperature approximation to obtain the external 
admittance of a cylindrical plasma capacitor. The appropriate differential 
equation was solved after using inner-outer expansions connected through 
the resonance region. A WKB description was used for the inner region, 
where resonant waves are essentially evanescent, while a travelling wave 
description was used for the outer region, where Landau damping is important. 
Because the wave nature of the solutions was assumed in the outer region, 


2 



this theory is appropriate only for higher order resonances. Similar 
difficulties have been shown to occur in the simpler one-dimensional model 
(Harker, Kino, and Eitelbach, 1968; Miura and Barston, 1971; Peratt and 
Kuehl, 1972). 

Variational methods offer an attractive alternative to these treatments. 
They have been used previously to estimate plasma resonance frequencies 
with simplified trial functions. Resonances of a cold inhomogeneous plasma 
were treated by Crawford and Kino (1963). Using the variational principle 
established by Sturrock (1958), Barston (1963) approximated the dispersion 
relations for wave propagation along an infinite cold plasma slab, and 
along the interface between two semi-infinite, counter-streaming cold 
plasmas. Some general features of the guided waves on a cold, transversely 
inhomogeneous plasma column in an axial magnetic field were studied by 
Briggs and Paik (1968). These papers (Crawford and Kino, 1963; Barston, 

1963; Briggs and Paik, 1968) show that, with appropriate variational 
principles and judicious choices of trial functions, useful results can 
be obtained with relative ease by the variational approach. 

A theoretical variational formulation for the electrostatic resonance 
oscillations of a warm, inhomogeneous plasma column in a dc electric or 
magnetic field of arbitrary direction was presented by Barston (1965) 
with the adiabatic index, y , taken as unity. The variational prin<- 
ciple to be presented here, however, is not restricted in the values of 
y . One important feature in Barston’s (1965) analysis is that the 
rf electric potential was treated as the solution of the rf Poisson 
equation, with the rf electron density considered given. It will be 
seen that the coordinate functions to be used here are defined in a 
similar fashion. A variational method of the Rayleigh-Ritz type has 
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been applied successfully by Dorman (1969) to a one-dimensional, warm, 

and field-free plasma with arbitrary dc density profile. A single 

second order differential equation for the electric field was obtained, 

and shown to have hermitian operators. The variational principle to 

be used here differs from Dorman's (1969) in that we are dealing 

directly with a system of Euler-Lagrange equations. In so doing, we 

can keep down the order of the equations, and are able to consider 

warm inhomogeneous plasmas in more than one dimension. 

In this paper, we shall show that by appropriate definitions of the 

rf equation of motion and pressure term, the equations of the hydrodynamic 

model used by Parker, et aJ. . (1964) become Euler-Lagrange equations. This 

will enable us to demonstrate the effectiveness of the Rayleigh-Ritz 

procedure in estimating the resonance frequencies of an inhomogeneous plasma. 

The associated numerical method mainly involves evaluations of definite 

integrals and solutions of finite algebraic eigenvalue equations, and is 

2 2 

applicable over the entire range of a /\ » 1 for estimating the first 

few resonance frequencies. 

In most of the papers that deal with the electrostatic resonance 

problem (Crawford and Kino, 1963; Parker et al., 1964; Harker et al., 1968; 

Baldwin, 1969; Dorman, 1969; Miura and Barston, 1971; Peratt and Kuehl , 

1972), it is assumed that the rf plasma current normal to the glass wall 

2 ,~~2 

is zero. However, in the low temperature limit, a -» co , the main 
resonance frequency seems to agree with that of cold plasma theory, in which 
the normal rf plasma current is retained. We consider this problem and 
show that the resulting difference in predicted resonance frequencies is 
negligibly small for low pressure positive columns. 
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In §2, we present the basic equations, the corresponding Lagrangian, 
and the procedure to be applied in the variational approach. In §3, the 
numerical methods are explained before comparing computations with those 
of Parker ej; al . (1964). The paper concludes with a brief discussion 
in §4. 
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2 . THEORY 


For a low pressure positive column, moment equations with scalar 
pressure and negligible heat conduction are appropriate when the wave 
phase velocity is much larger than the thermal speed. For the first few 

electrostatic resonances, the wave phase velocity may be scaled to $ a , 

2 ^ 
where U) p [ = (e n^O/me^ ' ] is the axial plasma frequency. Thus we require 

2 7~2 

a /)'- D » 1. A stationary ion background will be assumed, since we are 
interested only in electron resonances. Dissipation due to collisions, 
and Landau damping, will be neglected . Our analysis will consequently be 
valid only for the first few resonances. Also, the analysis will be quasi- 
static. Apart from some differences in definition of the rf equation 
of motion and pressure term, the equations we shall use are essentially 
those used by Parker, et al . (1964). The equations are generalized here 

to include dc- magnetic field, B , and electron drift velocity, v 

— (J ' -L0 

We have, 


at 


+ v * (nv) = 0 , 


€ 0 ?-£ + e ( n " n r ) = 0 , 

mn ("ft + + VP + en (£ + v; X B) = 0 jat r Q + . ( 1 ) 

Specialized to small perturbations, these reduce to the dc equations, 


= 0 • 


V 5 o + e ( W = 0 - 


mn oV 9 .V VP o + e VVVV = 0 < at V - 

and rf equations, 

Sn i 

— + <7-< Vj + “j.V = 0 . WZl + en 1= 0 , 

+ C^-£)vP 0 + VP X + £.?(v p ) 

+ en 0 ( i " 7 ~o + ii + i x * 0 + x V = 0 


( 2 ) 


(3) 
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where £ = d^/dt = B|/dt + . In these equations, m and -e are 

the electron mass and charge; n and P are the electron density and 
pressure; E is the electric field; is the ion density; is the 

vacuum permittivity; and ^ is the perturbation displacement for the 
electrons (Figure 1). 

The magnitude of is relatively small in the plasma region, but 

increases in the sheath region from the ion-acoustic speed to roughly the 
electron thermal speed at the glass wall (Self, 1963; Parker, 1963). We 
shall consequently neglect it in our analysis. However, due to the presence 
of non-zero v^, jn^, and an rf electric field at the wall, a non-zero rf 
normal current term arises, and hence an rf surface charge term. The 
inclusion of this surface charge term, in the cases where the electron 
rf excursion exceeds the Debye length, is equivalent to the use of the 
dielectric model for a cold plasma column (Crawford, 1965). Further 
discussion of this surface charge term will be given in 43.4. With v 
neglected, the following relations become appropriate 


Vt* = n o (0)f( £ ) • 


"i = 


f(r) = exp[-ecp Q (r)/KT e ] , 

II = 91/3* - 


(4) 


where k is the Boltzmann constant, cp (r) is the dc column potential, 

E^r) - -wp 0 (r) , (5) 

and the first equation of (3) has been used to obtain n^. The rf 
electron pressure, P^ , is determined by the adiabatic equation of state, 

p (£o + I>/ p o<£ 0 ) = n( £o + i ) 7 / n o ( £o )7 ■ p o = V T e ' (6) 

The form of (6) can be understood by reference to Figure 1, and 

follows from the fact that when a cell is displaced from r to r + £ 

-0 ~0 1 
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P is defined for the given cell t rather than for the given position, , 
The adiabatic equation of state must consequently be applied to the same 
cell, before and after the displacement. Using the usual definition of 
perturbations , 

= P o ( I ) + P !<£) . n(r) = n Q (r) + n^r) , (7) 

and (4), we obtain 

p i = '?V1- I' VP o • < 8 > 

The rf force law in (3) is obtained by comparing the force laws in (1) 
and (2) in the same fashion (Newcomb, 1962). 

It is now straightforward to use (4), (5), and (8) to rewrite (3) 
in terms of only £ and the rf potential, cp 1 (-^p 1 = E )., 

mn oL - ' v I' vP 0 ' ^ p 0 ,7<v ‘i ) 

- en 0 („p 1 + J-VVq> 0 - jxV = 0 , 

2 

e Q V cp 1 + ev(n^) = 0 . ( 9 ) 

Equations in (9) can be Fourier-transf ormed , normalized, and expressed 
in a cylindrical coordinate system, (r,0), for a column of cylindrical 
symmetry. 


fi 2 f£ r - n c nf? 


e 


A D [ 7 ^r + R^r R Se ) + R^e 1 


+ ^d ^ ( ^ + i^R C) ' + 1(i o^ +i, i ) = 0 ’ 


( 10 ) 


n 2 fC e - n c fi*e r + aV |(tc r -Ce) +■ t-l y* |(c; + i 


r - - r 

h r R 


f (i*o C e + 1 * 1 ) = 0 • 


<u) 
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( 12 ) 


F 1 1 + <“C P >' f £ e - 0 ■ 

where the derivative with respect to R is denoted by (') and 

V°> = 0 ’ i = ir £ r + ie <iC 0 ) ’ 

CO 

« 1 (R,0,T) = 

■t = ~CO 

with i denoting unit vectors. The normalized quantities are defined as 
R = r/a , T = m t , O - Oj/dD , 

Ir r 

n c = eB 0 /r ltUp , C = |/a , 

A D " > d ! &2 = £ 0 KT e/ n o (0)e2a2 ’ * = <pVV 0)ea2 > (14) 

and a static axial magnetic field, , has been included, 

2 , 1 Lagrangian Density 

The forms of the force law in (3), and the rf pressure in (7) , 
represent the important differences from the paper by Parker, _et aJ. (1964) 
in that they make (9), as well as (10)-(12), systems of Euler-Lagrange 
equations without having to restrict the values of y (Barston, 1965). 

In one of the models used by Dorman (1969), an appropriate pressure term 
similar to (7) was used without the benefit of the rf force law in (3). 

As a result, he was able to establish the variational principle only for 
the one-dimensional case. 

The Lagrangian corresponding to (10)- (12) can be shown to be the 
following, by straightforward application of Hamilton’s variational principle, 


/ 


(R)exp 


i(QT+te) 


(13 ) 


\ 


10 



*2 = 


a 

tf 


dQ£ 2 (Q f t) , 


Z f AQ.,l) = OA^QB+H , 
m c 


t=~ eo 0 
1 


A = / fRdR [\Q\ + 


(|C r | 2 * k 6 | 2 ) - B -f fRdR (Ve + £* r C 0 ) ■ 


H = ( 7 T 4 -T ' ) - V 0 S - I+F+F' < 0 , 


1 

T = j f RdR| £ ; + I c r - I g 2 ♦ I [fc*(c; + c r - ^c e > + C - C -] R= 1 ■ 

Tf 

1 

T' = J fRdR [i c;ac*-c*) + s CJ«-$ - C*e> + »•«>■] 


+ I [ f ^«Ve r > + f CS«W + °- c -] R= i ’ 


1 

S= J fRdR (*o|C r | 2 + | *okJ 2 ) ’ 

0 

1 1 

I = - f fRdR + £ C e J* + c.c.j , F = J RdR (|i ^| 2 + ^ |§ x | 2 

0 

R 


iv, ** 

b t C 

F ' =y € g RdR (|$'| 2 + jij 2 ) + J RdR (|$'| 2 + ^ |*J 2 ) , 


R b = - , r = £ 

b a c a 


( 15 ) 
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where c.c. denotes complex conjugate. Equation (15) is appropriate 

to the configuration shown in Figure 2 of a concentric metal cylinder 

surrounding a glass tube, of relative permittivity e , that contains 

S 

the plasma column. In the expressions for T and t' , the boundary 
terms are included to modify the natural boundary conditions on £ 
an d Qq (see for example, Courant and Hilbert, 1953) that would other- 
with be unphysical. In the expressions for H and S , 0^(r) anci V ^ 
are defined as 

V r) = - V r) / v o - v 0 = - V a) • 

Substitution of the exact solutions of (10)- (12) would make £„(£!,£) 

zero. We see from (15) that, for negligible Q , A , and V , the 

c D 0 

resonance frequencies are determined essentially by the values of I, F, 
and F / . Since is approximately proportional to (Self, 1963; 

Parker, 1963), the effect of higher electron temperature is to raise each 
resonance frequency. When ft ^ 0 , the roots of <2, (0,-t) = 0 are 

C 

2 = fl c B/2A ± [(Q c B/2A) 2 - H/A . (16) 

Since Q and - Q are indistinguishable in experimental observations, we 
see that all the resonance frequencies are predicted by (16) to split 
in two, in agreement with the theoretical results of Barston (1965), and 
Vandenplas and Messiaen (1965). For sufficiently small Q , when the 
values of B and A are not greatly affected by the presence of an 
axial static magnetic field, the amount of the split, |q B/Aj , will 
be proportional to Q . This is in agreement with the observed splitting 
character of the main dipole resonance frequencies (Crawford, Kino, and 
Cannara, 1963; Messiaen and Vandenplas, 1962). Furthermore, since 
A 5: |B| , where the equal sign applies when £ = , this split is predicted 

to be always less than or equal to Q 
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2-2 Rayleigh-Ritz Procedure and Coordinate Functions 


For a single linear Euler-Lagrangian equation, the Rayleigh-Ritz 
procedure is efficient in obtaining an approximate solution, by use of a 
weighted summation of a set of judiciously chosen coordinate functions. 
These coordinate functions must be linearly independent and complete, 
and satisfy the boundary conditions specified by the problem. The 
weighting coefficients that appear in the approximated Lagrangian are 
varied independently. This results in a system of algebraic equations 
that take the place of the original equation. Theoretically, better 
approximations can be obtained by using more coordinate functions. When 
eigenvalues are involved, the approximate eigenvalues always converge to 
the exact values from above (see for example, Mikhlin, 1964). 

Much less attention has been given to the analogous problem for a 
system of linear Euler-Lagrangian equations, e.g. (10)-(12), though the 
theoretical extension of the variational method to a system of second order 
elliptic equations has been mentioned by Mikhlin (1965); the way to set 
up the corresponding coordinate functions is similar to that for the single 
equation case, i.e. the coefficients of each dependent variable are varied 
independently . 

For our problem, in which (10)- (12) are not in elliptic form, the 
coordinate functions and the coefficients must be more restricted. In the 
Appendix to the paper it is shown that acceptable estimates of resonance 
frequencies and eigenfunctions can be obtained for our problem provided 
that the coordinate functions chosen for each dependent variable are 
related by (10)- (12), 

By expanding Q , Q , and $ in power series of R, and substituting 
r y l 

in (10)- (12), we see that for small R , oc ^ cc R „ f and $ cc R^' , 
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Thus, for JL 2 . 1 , the solutions are well behaved at R = 0 . If we 
also choose even functions for f(R) and §^(R) , then and are 

even in R , and $ is odd in R , when t is odd, and vice versa . 

Since there are no other singularities in (10)-(12), polynomials in R 
constitute appropriate coordinate functions for our problem. For convenience, 
the coordinate functions chosen for r will be 

tor 

r . = R^ -1 - R^- 1 (j = 1,2,...) , <17) 

^rj 

which conform to the usual assumption of zero normal rf current, since 


C -CD = 0 . 

Rather than choosing r and . $ , independently, we must determine 

them via the original differential equations and (17). After eliminating 
§ in (10) and (11), it follows by using the second equation of (4) that 


+ ^ = t R n c Q + WrUf o ]( "C + C r '-£-6 9 > ■ a«) 


With a given expression for $^(r) , and an assigned value of ft , e.g. 

H=l, Q can then be easily determined for any given r 

wj ■ rj 

An immediate question arises concerning the dependence of the resulting 

variational estimate of resonance frequencies on the size of Q chosen 

arbitrarily here. We have found that the first few resonance frequencies 

-4 

do not change by more than 10 u> when ft in (18) changes from 0.4 to 1 

for all the values of l/Ajj used in this paper. If this were not the 
case, an iterative procedure would have to be used, i.e. the resulting 
variational estimate of ft would have to be used in (18) to obtain a 

new set of . These would be used in turn to obtain improved frequency 

estimates, and so on. 

The corresponding coordinate function, $ . , is obtained by solving 

i , J 

(12) (Barston, 1965), with conditions of continuity of potential and normal 
displacement across the boundaries defined in Figure 2. According to the 
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discussion given in the Appendix, identical coefficients are assigned to 
each set of coordinate functions, 


N N N 

=E ’ C 9 = £ a.£ ej , • («> 

j=l j=l 3=1 

before substitution in the Lagrangian, £ (0,1), of (15). The resulting 
Lagrangian then gives the algebraic Euler-Lagrange equation below, 



(ft 2 A ,, -’ffiB + H ) a . = 0 
Ji c ji ji' j 

1 

fR dR(£ ri C r j + C-eiCej) » 



(i = 1,2, . . . , N) , 


( 20 ) 


where A. . , B. . , and H. . are the matrix elements of the integrals, A, 

B, and H, given in (15), respectively, and are obtained by substituting 

the coordinate functions in a fashion given by the above A exoression 

ij 

Equation (20) can be transformed into a generalized eigenvalue problem 
(Dorman, 1969) 

2N 

- 

✓ , (^A " B )b = 0 (i = 1,2,... , 2N) , 


b . = (a. >Qa . ) 
J i l 


U = 1,2,..., 2N, i = 1,2,... 


N) , 



with superscript + signifying the transposition of a matrix. Equation (21) 
is now solvable by standard computer codes. 
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3. NUMERICAL METHODS AND RESULTS 

The computing procedure is straightforward: 

(i) read in physical and computational parameters, 

(ii) compute coordinate functions .. , and all other functions 

appearing in the Lagrangian , » at intervals Ar , 

(iii) compute (ft,£) by Simpson's rule to obtain A. , , B . . , and 

Z Ij lj 

H . . , and 

(iv) solve (21) for f2 and a_. , and compute relevant eigenfunctions for 
N, N-l, and N-2 coordinate functions. 

3.1. Approximate DC Density Profile 

A density profile which approximates Parker’s results (Parker, 
1963), and is convenient for both analytical and numerical manipulation, 
is given by 

f(R) - exp [- \4> 0 (R)J , <|> 0 <R) = gR 2 + (l-fS)R h , (22) 

2 

where h(> 2) is an even integer, S < 1 , and T\ [= V /A ] is the wall 

'w O' D 

value of the potential function, 7|(r) , used by Self (1963) and Parker 

(1963). The particular form of (22) is used to justify the choice of 

coordinate functions of (17). Furthermore, with the use of (22), 

and § can now be solved analytically in terms of power series in R, 

* j J 

in addition to the numerical solutions of (12) and (18). Comparison of 

the solutions by the two methods will provide an estimate of the degree 

of accuracy achieved in obtaining £„ . and $ 

H J 1 , J 

The values of p and h are varied until f(r) best approximates, 

by least-square deviation, the profile given by Parker for specified values 

of Aq and T] . The resulting profiles are shown in Figure 3, and the 

corresponding values of j3 and h are given in Table 1. It will be 
, 2 

seen that as 1/A^ increases in value, (22) decreases in accuracy because 
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FIG. 3. Comparison between (22) (indicated by dots), and Parker’s 
density profile, for the parameters given in Table 1. 


Table 1 


Parameters used in (22), and Figure 3, for a mercury- vapor 

2 2 

positive column. The conversion between l/A and 1A Q 
is obtained by the calculations of Parker (1963). 



VA o 

VA D 

T] 

‘w 

3 

h 

* 2 
SA 

+ EA 2 R 

a 

3.4 x 10 2 

7.2 X 10 1 

6.72 

0.447 

4 

* 

6.2 

x io' 4 

b 

1.3 x 10 3 

5.1 x 10 2 

6.60 

0.221 

8 

* 

1.3 

x io 4 

c 

8.2 x 10 3 

4.5 x 10 3 

6.52 

0.159 

20 

1 6 . 2 

-4 

X io 

d 

6.7 x 10 4 

4.3 x 10 4 

6.44 

0.142 

54 

t2 .9 

X 10“ 3 

e 

6.2 X 10 5 

4.2 x 10 5 

6.40 

0.141 

156 

f7 . 5 

X 10' 3 

f 

CD 

00 

1.08 

0.644 

12 

1 1 . 3 

X lo' 4 


L9 



of the increasing degree of steepness displayed by the density profile 

in the sheath region. At zero temperature, the sheath is omitted and 

an accurate approximation can again be obtained . The approximation of 

(22), however, was found to be sufficient to give accurate frequency 

2 

and eigenfunction estimates for all of the values of 1 /a listed in 
Table 1. 
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3.2. Solutions for Q and $ 

^ J > J 

By the use of (14) , (15), and (19) , can be put in the 

form of a rapidly converging series, 


- R 


i- 


1 - E vm 


ft+2j-t(i-l) (h-2) -1 


i=l 


Yl " -2 




ft^+(&+2j-2)|-S- - 2$ (7-1) V, 


ft 


j- (JU-2J-2) + SMI-2 & (7 -IHVq 


<7"l)(l~P)hV 0 [^a J 1 -a+2j-2)] 

^ j , 2 = ^2 1 

^2 (t+2 j+h-4)+Q c Q-2p ( 7 -l RV q 
'C 


S, Z (^l)(l-0)hV o 


J*1 


ifi 2 (h-2) 




(i > 3). 


(23) 


Since the maximum value of Vq of interest is roughly 0.02 and 

ft 't 1. , one needs at most five or six terms in (23) to attain a precision 

of 10“ 8 for C a> . 

6j 

To solve for <l> , , the predictor-corrector method of Adams-Bashf orth 

1 » J 

(Fox, 1962) has been used on (12), which can be reduced to the form, 


y l = y 2^ R > y 2 = 8^ + i2y l^ R • 


y l = *1,3 ’ 


y 2 - R *l,j • 


g(R) = f(-RT] w $^ rj + R C ;. + C r . - IQ ) . (24) 
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N N L 

Since the complementary solution of (24) is $ . = y = R , and its 

I > J 

P , A ^ 

particular solution, y , is proportional to R near r = 0 , 

p > p 1 

y (0) and y (0) are both equal to zero. So the starting values of 
1 2 

P P 

y and y are well-behaved, and easily obtained by Taylor series 
1 2 


expansion near R = 0 . The total solution of $ 

y J 


can then be written 


as 


*i.j -°/ + ■ 


(25) 


where c . is determined by imposing the boundary conditions of . 

By making the interval Ar = 0.01, and using double precision, § . 

— 8 

can be calculated to within 10 « This is arrived at, first, by com- 

paring results that use different values of A r j anc ^ secondly, checking 
against solutions of (12) obtained by power-series expansions in R . 

3 . 3 Numerical instability 

In the process of solving the algebraic equation, (21), the size 
of N is limited by the inaccuracy involved in obtaining A , etc, 

•*- J 

This inaccuracy introduces a numerical instability whenever the coordinate 
functions are not orthogonal functions with respect to the differential 
operators of (10)-(12) (Mikhlin, 1971, Chap. 2). The situation is best 
illustrated by an example in which the dipole resonance frequencies 
corresponding to Figure 3(b) are calculated for different values of N 
while holding the size of A r constant at 0.05. 

As shown in Table 2, as N is increased from 2 , the first few 
resonance frequencies are approached from above with rapidly stabilized 
estimates, When N is increased beyond 8, undesirable fluctuations 
larger than 10 4 , and clearly erratic changes in the values of Q , start 
to appear. In the case N = 9, for example, one would obtain an erroneous 
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Table 2 


Dipole resonance frequency estimates obtained ^ith A r = 0,05 
for the case of Figure 3(b). Significant numerical insta- 
bility sets in when N a 8 . The ‘best estimates’ are 
obtained with A r = 0.01 and N = 9 . 


N 

r 

Main 

First 

Second 

Third 

Fourth 

Best estimates: 

0.4447 

0.7099 

0.8746 

0.9924 

1.080 

2 

* " 1 

0.4703 

0.8813 




3 

0.4486 

0.7446 

0.9250 



i 

4 

0.4450 

0.7155 

0.8818 

1.022 


\ 5 

1 

0.4448 

0.7119 

0.8757 

0.9958 

1.084 

j 6 

0.4448 

0.7117 

0.8753 

0.9956 

1.084 

i 

' 7 

0.4448 

0.7104 

0.8751 

0.9947 

1.082 

1 8 

i 

0.4448 

0.7098 

0.8748 

0.9921 

1.081 


\ 
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fundamental resonance frequency. Characteristic of the variational 
nature of the Lagrangian, (0,-6), more serious errors are found in 
the approximate eigenfunctions, than in the resonance frequencies. The 
optimal combination of N and &r , that produces acceptable results 
in the shortest computation time, can be obtained by trial and error. 
Repeated solution of (21) for a few adjacent values of N thus becomes 


an economical technique: this requires computation of the matrices A. . , 

t J 


etc. only once, and offers safeguards against obtaining erroneous results 


due to numerical ins tabilities . 


3.4 Computer Results 

As a practical example, the approximate density profiles given 
in Table 1 have been used to predict dipole resonances for Tube No. 1 


used by' Parker, et al. (1964) (£, = 1; a = 0,5 cm; effective relative 

permittivity at the surface of the column K ^ =■ 2.1). The 

2 

computation time varies roughly as N . With N = 10, a typical calcu- 
lation takes about 40 seconds, and requires a core space of less than 
100K bytes in an IBM 370/67 machine. 

The resulting approximate solutions for ? , i£ , n . co , and cd / 

®r 1 Y 1 T 1 

are plotted in Figures 4-6. The density and the radial electric field 
solutions, n^(R) and ep^ resemble very closely those given by Parker, et al . 
(1964), and Parbhakar and Gregory (1971), respectively. The relative 


amplitudes of | and §_ shown in these figures are retained, revealing 

I* Tj 

2 

that as l/A increases, | ^ progressively dominates over | . For 

‘ D r 0 

~~2 

1 /([ > 4500 , it will be seen that the perturbations should be progressively 

2 

compressed toward the sheath region as l/A increases in value. These 
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solutions are not reproduced here because they also exhibit undesirable 
oscillations with wavenumber equal to N , an expected characteristic 
when we try to approximate rapidly varying functions with truncated 
polynomials. Since only a moderate computer storage is used for N £ 10 , 
there is room to increase N, and decrease , to obtain better approxi- 
mate solutions. However, this is considered unimportant for our purpose, 
since we are able to obtain good frequency estimates for this region with 
N s 10, as Figure 7 reveals. 

The corresponding estimates of resonance frequencies are shown in 
Figure 7. Since the electron temperature corresponding to the experimental 
resonance data of Parker et al. (1964) was adjusted to fit their theoretical 
spectrum, it would be reasonable for us to make a similar adjustment. As 
is evident from Figure 7, however, no such adjustment is necessary. Indeed, 
our result seems to be in slightly better agreement with the T g = 3 eV 
data. The minor differences between the two theoretical results probably 
come from the differences in the rf equation of motion and pressure 
term used in the two treatments 

Similar to other papers (Crawford, 1964; Parker et al. 1964; Harker 
et al. 1968; Baldwin, 1969; Dorman, 1969; Miura and Barston, 1971; Peratt 
and Kuehl, 1972), we have assumed zero normal rf plasma current density 
at the glass wall, through the form of £ in (17). This is appropriate 
when the plasma is sufficiently warm that the electron excursion velocity 
is much smaller than the thermal speed, and f(l) « 1 . This assumption, 
however, is inconsistent with the dielectric model for a cold plasma 
column, where normal rf plasma current must be included. It is of 
interest to ask why the main resonance frequency of a warm plasma column, 
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in the limit of low temperature, approaches the principal resonance of 
a cold plasma column. To answer this question, we need only change (17) 
to 

C rJ = tf 1 " 1 + R t+2J_1 (j = 1,2,...) , (26) 

and impose the requirement of continuity of normal displacement in the 
form 

$'<!“) + f(DC (1) = e $,'<i + ) . (27) 

i *r g i 

The resulting solutions of are found to be only slightly different 

from the previous case near R = 1 [Figures 4-6, where the dashed lines 

correspond to the use of (17)]. Furthermore, the main resonance is 

~2 

lowered by less than 1 per cent for all of the values of 1/A^ used 
here, including the case l/A^ -* . This is well within the experi- 

mental errors. 


\ 
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4. DISCUSSION 


In this paper, we have applied thq Rayleigh-Ritz procedure to a system 

of three Euler-Lagrangian equations that describe the electron resonances 

of a nonuniform warm plasma column. It is shown that accurate frequencies 

for the first few resonances can be obtained for the entire range of 
o 

l/^D ' > ' > ^ • Results which agree closely with those of Parker, ejt al . 

(1964) have been obtained 

Contrary to the case of a system of elliptic equations, where the 
coefficients are assigned independently to each dependent variable 
(Mikhlin, 1971), we have found that for (10)-(12), the same coefficient 
must be assigned to each set of coordinate functions, e.g. (19). In 
addition to the usual requirements , that the coordinate functions must be 
linearly independent and complete, we have chosen that they be set up in 
accordance with (10)-(12). 

The present method can be easily modified to include the effects 
of electron dc drift, dc magnetic field, and ion motion. With the axial 
dimensions and rf magnetic field included, this procedure would be efficient 
in solving travelling wave problems in a nonuniform plasma waveguide. 
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APPENDIX 


Here we shall show that the coordinate functions, , and 

4 . , must be assigned the same coefficients, a , for the Rayleigh- 

1 , J J 

Ritz procedure applied to £ (Q,£) of (15) to be successful: it will 

be shown that, by restricting these coordinate functions according to 
(10)- (12 ) , the appropriate Rayleigh-Ri tz procedure can be established. 
A system of second order elliptic equations can be written as 


( V u k> ' + AE jk u ■ 0 > V - D kj * 


E 41 = E , 

jk kj 


(A.l) 


th ' 

where U, is the k dependent variable, the matrices D and E 

k J K J K 

are functions of the independent variable x , the eigenvalue is X(^ 0) 5 
and the summation convention has been used . Ellipticity demands that 


=*• 0 > E jk a j“k 2 0 ■ 


(A. 2) 


for any real non-zero vector o;^ . The solutions of (A.l) then admit 
of variational estimates, as outlined by Mikhlin (1965). 

When we apply the Raylei gh-Ri tz procedure to the Lagrangian for (A.l), 


■/ 


dx (D^U'U,' - X E . . U . U. ) , 


jk j k 


J* j ^ 


(A. 3) 


the coefficients preceding the coordinate function for each U can be 

varied independently. Suppose a set of legitimate trial functions, 

n = C V. (summation convention not used here) with coefficient C, , are 
JK K rL f* 

used in (A. 3). Because of (A. 2), we can obtain crude estimates of \ 

larger than the lowest eigenvalue, even if only one of the C is non-zero. 

If the same procedure is applied to the Lagrangian in (15), we will 

obtain erroneous estimates of Q . Consider the case of a cold, uniform 

plasma column with = 0 , so that H=-I+F+F / . Suppose C , 

0 r 
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C , and C are the coefficients for the coordinate functions £ . , 
0 0 

r , and ®. . . respectively. Then if C = 0 , the value of the Q 

b r0 i, J 0 

estimate will be zero, and fail the requirement of the Rayleigh-Ritz 


,V2 


sequence . 

Consequently, to make (-H/A)" 7 " from (15) at least non-zero, we 

must use a single coefficient for each set of the coordinate functions 

r . , r . , and $ . . Merely making C = CJ = = a is insuffi- 

cient to produce a legitimate Rayleigh-Ritz sequence, because the value 

of $ , for example, can be arbitrarily small in comparison with 

1 > J 

r and r , making the Q estimates also arbitrarily small. We 
b rj ^Oj 

see that restricting these coordinate functions according to (10)- (12), 
is sufficient to reduce the resulting approximated Lagrangian to a single 
Euler-Lagrange equation. The appropriateness of the resulting Rayleigh- 
Ritz procedure can then be guaranteed. It is conceivable, of course, 
that there may be less restrictive choices of appropriate coordinate 
functions corresponding to the Lagrangian in (15), but we have not chosen 
to pursue this point. 
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